function [X, res] = monmerplt(y, handles, x1, x2)
% subfunction to display calculated and experimental data
% Y        calculated saturation recovery data
% handles  structure with handles and user data
global hstart 
global dat
global deltao

% calculate x-axis in micro-seconds
X(x2-x1+1) = 0.0;
for ix = 1:(x2-x1+1)
  X(ix) = (hstart + (ix-1)*deltao)*1E6;
end
xlabel('Time (microseconds)');

% initialize arrays
% res = residuals
% pct = scaled calculated data
% pct2 = scaled experimental data
res(x2-x1+1) = 0.0;
pct(x2-x1+1) = 0.0;
pct2(x2-x1+1) = 0.0;

% The data is scaled as a percentage of 110% of its range.  Using
% 110% allows the maxium value to be 10% below the top gridline.
scale = 1.1;
Ymax = max(y);
Ymin = min(y);
Ymid = (Ymin + Ymax)*0.5;
Yrng = (Ymax - Ymin)*scale;

% The maximum value of the experimental data is calculated as the
% average of the last 9 points.  The data is then scaled as a
% percentage of 110% of its range.
Y2max = 0;
for ix = x2-8:x2
  Y2max = Y2max + dat(ix);
end
Y2max = Y2max/10;
Y2min = min(dat);
Y2mid = (Y2min + Y2max)*0.5;
Y2rng = (Y2max - Y2min)*scale;

ssq = 0;  % sum of squares of residuals

for ix = 1:(x2-x1+1)
  pct(ix) = (y(ix)-Ymin)/Yrng;
  pct2(ix) = (dat(ix)-Y2min)/Y2rng;
  res(ix) = (dat(ix)-Y2mid)/Y2rng - (y(ix)-Ymid)/Yrng;
  ssq = ssq + res(ix)*res(ix);
end
disp(['ssq: ' num2str(ssq)]);
set(handles.ssqNum, 'String', ssq);

hold off
xlim([x1 x2]);
plotyy(X, pct2, X, pct);

axes(handles.residPlot);
plot(res);
end
